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Abstract 

Impedance matrices are obtained for radially inhomogeneous structures using 
the Stroh-like system of six first order differential equations for the time harmonic 
displacement-traction 6- vector. Particular attention is paid to the newly identified 
solid-cylinder impedance matrix Z(r) appropriate to cylinders with material at r = 
0, and its limiting value at that point, the solid-cylinder impedance matrix Zq. 
We show that Zq is a fundamental material property depending only on the elastic 
moduli and the azimuthal order n, that Z(r) is Hermitian and Zq is negative semi- 
definite. Explicit solutions for Zq are presented for monoclinic and higher material 
symmetry, and the special cases of n = and 1 are treated in detail. Two methods 
are proposed for finding Z(r), one based on the Frobenius series solution and the 
other using a differential Riccati equation with Zq as initial value. The radiation 
impedance matrix is defined and shown to be non-Hermitian. These impedance 
matrices enable concise and efficient formulations of dispersion equations for wave 
guides, and solutions of scattering and related wave problems in cylinders. 

Introduction 

Impedance provides a useful tool for solving dynamic problems in acoustics and elasticity. 
A single scalar impedance is usually sufficient in acoustics, whereas a matrix of impedance 
elements is required to handle the vector nature of elastic wave motion, particularly in the 
presence of anisotropy. The use of impedance matrices can offer new insight because their 
properties are intimately related to the fundamental physics of the problem, as, for in- 
stance, the Hermitian property of the impedance matrix which is directly linked to energy 
considerations. A classical example is Lothe & Barnett's surface impedance matrix [l|,l2[ 
which proved to be crucial for understanding surface waves in anisotropic homogeneous 
half-spaces, with the result that it provides perhaps the simplest method for finding the 
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Rayleigh wave speed. Biryukov js], 0] has developed a general impedance approach for 
surface waves in inhomogeneous half-spaces based on the differential Riccati equation, see 
also (11 . Direct use of the impedance rather than the full displacement-traction wave field 
provides an efficient and stable procedure for computing high-frequency dispersion spec- 
tra. Several numerical schemes for guided waves and scattering in multilayered structures 
have been developed on this basis, e.g. jl, 0, |8|. These involve the 3x3 impedance matrix 
often called the 'surface impedance', although it actually differs from the 3x3 surface 
impedance of a half-space. It is useful to further distinguish the familiar 3x3 ('condi- 
tional') impedance from a 6x6 ('two-point') matrix more closely related to the matricant 
of the system equations. The nature of these impedances have been analyzed using the 
Stroh framework for homogeneous and functionally graded plates jH, [lO| . It is noteworthy 
that both impedances are Hermitian under appropriate physical assumptions; however, 
their Hermiticity implies a somewhat different energy-flux property than the Hermiticity 
of the half-space impedance. A bibliography on the impedance matrices for piezoelectric 
media may be found in 0, [HI . 

The above review concerns rectangularly anisotropic materials and planar structures. 
The objective of this paper is to provide an equally comprehensive impedance formalism 
for time-harmonic modes of nth azimuthal order in radially inhomogeneous cylindrically 
anisotropic materials of infinite axial extent and various circular configurations. An im- 
portant element in this task is the Stroh-like state-vector formalism developed for such 



materials in [12|. The results of [12[, which are based on the matricant in a Peano- 
series form (particularly the definition of the 'two-point' impedances similar to the case 
of planar structures) are however only relevant to a cylindrical annulus with no material 
around the central point r = 0. The intrinsic singularity of elastodynamic solutions at 
the origin of the cylindrical coordinate system, which rules out the Peano series, is an 
essential distinguishing feature as compared to the Cartesian setup. The problem can be 
readily handled in (transversely) isotropic homogeneous media with explicit Bessel solu- 
tions; however, it becomes considerably more intricate for cylindrically anisotropic and 
for radially inhomogeneous solid cylinders. The main analytical tool in this case is the 
Frobenius series solution. The milestone results on its application to hornogeneous and 
layered cylinders of various classes of cylindrical anisotropy include [H, [13, [l5 , [H, [13, [H 



(see also the review p^). The state- vector formalism based on the Frobenius solution 
for the general case of unrestricted cylindrical anisotropy and arbitrary radial variation 
of material properties [10] is of crucial importance to the present study. Another vital 
ingredient is the differential matrix Riccati equation for an impedance [3]. To the best 
of the authors' knowledge, this equation has only recently been used for the first time in 
elasticity of cylinders by Destrade et al. 2l|| who numerically solved it for an elastostatic 
problem in tubes. 

The presence of the special point r = distinguishes the solid-cylinder case from its 
Cartesian counterpart in many ways. Apart from the usual radiation condition at infinity, 
a similar kind of condition has to be applied at r = 0. The Riccati equation simultaneously 
determines the central-impedance at r = in a consistent manner while requiring it as the 
initial value for obtaining the solid-cylinder impedance. No other auxiliary (boundary) 
condition applies at r = 0, which is (again) unlike the surface or conditional impedance 
for, say, a traction-free plane y = 0. These observations point to the fundamental role 
of the impedance formalism in cylindrically anisotropic elastodynamics and actually call 
for a new type of the impedance matrix appropriate for solid cylinders. The concept, 
properties and calculation of the solid-cylinder impedance are among the main results of 
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this paper. 

The outhne is as follows. Background material on the matricant, impedance matrices 
and Riccati equations is presented in in a general context not specific to cylindrical 
configurations. In [|2] the governing equations for cylindrically anisotropic elastic solids 
are reviewed and the first order differential system for the displacement-traction vector is 
described. Some examples of the use of impedance matrices are discussed in ^ and in the 
process the solid cylinder and the radiation impedance matrices are introduced. Methods 
for determining the solid cylinder impedance are developed in ^ This section provides 
a detailed description of the Frobenius solution and its properties, and also discusses the 
Riccati solution. Both methods involve the crucial central impedance matrix, to which 
||5] is devoted, where explicit solutions are presented and general attributes delineated, 
including the important Hermitian property. The radiation impedance matrix is analyzed 
in ^ Explicit examples are presented in ^ for the central impedance matrix in different 
types of anisotropy, and the solid cylinder impedance is explicitly presented for transverse 
isotropy. Numerical results illustrating the Riccati equation solution method are also 
discussed. Concluding remarks can be found in ^ 



1 The matricant, impedance matrices and Riccati 
equations 

For the moment the development is independent of the physical dimension and the un- 
derlying coordinates. Consider a system of 2m linear ordinary differential equations 

The m— dimensional vectors U, V and the mxm submatrices Qj, j = 1, 2, 3, 4 all possess 
uni-dimensional spatial dependence on y, which may be a Cartesian or radial coordinate. 
The system matrix Q displays an important algebraic symmetry which is a consequence of 
a general flux continuity condition. The derivative of the scalar quantity ?7^T97, where su- 
perscript '+' means the adjoint (complex conjugate transpose) and T has block structure 
with zero submatrices on the diagonal and off-diagonal mxm identity matrices, can be 
identifled with the divergence of the flux vector P (to be deflned more speciflcally later). 
Thus, ■^(rj'^Tr)^ ~ divP, and hence, ([1]) implies the connection between flux continuity 
and symmetry of the system matrix [l0| 

Q = -TQ+T ^ divP = 0. (2) 

The vanishing of divP = assumes certain physical restrictions which will be described 
when the elasticity problem is considered in Section [21 

The 2mx2m matricant M(?/, y^) is a function of two coordinates deflned as the solution 
of the initial value problem 

^(y, Vo) = Q(l/)M(?/, yo), M{yo, yo) = 1(2^)- (3) 

The matricant may be represented formally as a Volterra or multiplicative integral eval- 
uated by means of a Peano series [22|, alternatively it may be expanded in a Frobenius 
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series [23|. Let rfa{y) {a = 1,2, ... , 2m) be a set of partial solutions, that is, a complete set 
of independent solutions of the homogeneous system ([T]), then M(?/,yo) = Mi!j)M ' (//o), 
where jV is the integral matrix (a first-rank tensor) jV (y) = ||?7i, ...,r/2m||- The propagator 
nature of the matricant is apparent from the property M(i/, yi)NL{yi,yo) = M(?/, yo), and 
in particular M(|/,|/o) = M(yo,y)~^, while the symmetry implies 

M(y,yo) = TM+(2/o,y)T. (4) 

Hence 

M-i(y,yo) = TM+(y,|/o)T, (5) 

that is, M is T-unitary [2^ . 

In solving problems one is often not interested in the individual fields V{y) and V(?/), 
but rather in their relationship to one another, and perhaps only at one or two positions 
such as boundary values of y. Accordingly we introduce the mxm conditional impedance 
matrix z defined such that 

V{y) = -tziy)V{y). (6) 

The conditional nature of this impedance arises from an auxiliary condition at another 
coordinate yQ [13], and may be understood from an equivalent definition of the matri- 

V M )-[m, m J (v (,„) j • "•>^^<= "^(f • = M J ■ 

Now suppose z(?/o) is the conditional impedance at y^, then 

V{y) = {M^-^M2z{yo))V{yo), 
V(|/) = (M3-zM4z(|/o))U(2/o), 

and the conditional impedance at y is therefore 

z{y) = t{M3 - 2M4z(yo)) (Ml - tM2z{yo)y\ (8) 

In practice, z(|/o) is often associated with boundary conditions on the level surface y = yQ. 
For instance, 'zero traction' and 'rigid boundary' conditions are specified by vanishing V 
and U, respectively, with conditional impedances 

^ fzMsM^i zero traction (V (yo) = 0), 
|iM4M2^ rigid boundary (U(yo) = 0), 

where M^- = {y, yo) in eqs. d?])-®. 

While it is possible to define the conditional impedance in terms of solutions of the 
2m X 2m linear system ([T]), the same system leads through a process of elimination to a 
uadratically nonlinear equation for the mxm matrix z: the differential Riccati equation 

dz 

- — h zQi — Q4Z — zzQ2Z — iQa = 0. (10) 
dy 

In this context the auxiliary impedance z{yo) serves as an initial condition a.ty = yo which 
once specified uniquely determines z{y) at other positions. The symmetry ([2])i renders 
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eq. flTU]) self-adjoint in the sense that if z is a solution then so is z+, which does not imply 
their equality. It does however imply that the differential Riccati equation flTII]) produces 
an Hermitian impedance, z(r) = z"'"(r), as long as the initial condition is Hermitian, 
z{yo) = z+(yo)- We will also find useful the algebraic Riccati equation associated with eq. 

(HDD, 

zQi — Q4Z — izQ2Z — iQs = 0, (11) 

the solution of which determines limiting values of the impedance, e.g. as \y\ 00, and 
can serve as the initial value for the differential equation ( !T0|) . 

We also introduce a 'two-point' impedance Z (y, yo) distinguished from the conditional 
impedance by its explicit dependence upon two arguments, and defined such that it relates 
the constituent parts of the 2m- vector at y and yo according to [1, [lO[ 



-V(y) ) - -^'^'V Ufe) 

Comparing ([6]) and ([12]) one might be tempted to surmise that the two-point impedance 
is composed simply of block diagonal elements Zi(?/,?/o) and Z4^{y,yo) identified as z(?/o) 
and —z{y), respectively, where z is the conditional impedance, and with zero off-diagonal 
blocks (Z3 and Z2). But the two-point impedance is more fundamental and thereby richer, 
as one can see by comparing ([7]) and ( !T2|) . implying 




M4M^^Mi-M3 -M4M^7 detMs 



(13) 



/Ml 


M2] 




CO 


mJ 


H 



— Z2 ^ 

ZZ3 — t'Zj^'Zi2 ^Z]^ — Z4Z2 



1 . . det Z3 

det M = -. 

det Z2 



where Mj = 'Wlj{y,yQ), Zj = Zj{y,yQ). The identity then implies the important 
properties that the two-point impedance is Hermitian, and that the matricant determinant 
is of unit magnitude, i.e., 

Eq. (ED^ Z = Z+, detM = e*'^ where = arg det(MiM4). (14) 

It follows directly from (|3])^ and Jacobi's formula that the phase satisfies the differen- 
tial equation ^ = — itrQ with initial condition (j){yo) = 0. The matricant is therefore 

unimodular (detM = 1) if trQ vanishes. Further properties of the impedance may be 
deduced by swapping the 'running' and 'reference' points y and yo in f|T2|) (i.e. in both 
U, V and Z), implying the reciprocal form 

-V(yo) ) - '^^^°'^Hu(yo) 
whence follows an obvious relation 

Z{yo,y) = -TZ{y,yo)T. (15) 
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The two-point impedance therefore has the structure 



Z 

7(„„\-l ^ M with Zi(?/,?/o) = -24(1/0, y), 

Z{y,yo)-\ I, with z^(^^^^)^_z3(^^^^). (16) 

As an alternative to eq. ([HD the conditional impedance at y may be expressed in terms of 
the impedance at yo by using the two-point impedance, 

z(y) = -Z4-Z3(z(yo)-Zi)^'Z2, (17) 

where Zj = Zj (y,yo)- Note that z{y) is Hermitian if z(yo) is. 

In the same way that the matricant M(t/, yo) satisfies an ordinary differential equation 
in y, viz. (|3]), it is possible to express the dependence of Z{y, yo) on y in differential form. 
Differentiating f lT2|) with respect to y and using ([1]) to eliminate the traction vectors yields 
an equation for the two-point impedance, 

^ + ZJi - J4Z + ^ZJ2Z + a3 = 0, where J. = ( ° ^ ) . (18) 

dy \y) ^j{y) ) 

The self-adjoint property of this equation is obvious because the two-point impedance is 
itself self-adjoint (Hermitian). Direct integration of the differential system (fT8|) subject 
to initial conditions at ?/ = yo is problematic because of the fact that all submatrices of 
Z(i/,yo) are of the form ±^(/J^ dy Q2)~^ as |y — yol — ^ 0, and hence undefined. Differential 
equations with well defined (finite) initial value conditions can be obtained for the block 
matrices Z~^ (j = 1,2,3,4) by simple manipulation of eq. (fT8l) . but we do not discuss 
this further here. It is interesting to note, however, that inspection of the block structure 
of (fT8|) shows that the equation for Z4 decouples from the other submatrices and it is the 
same as the differential Riccati equation ( ITOi) for the conditional impedance (under the 
interchange Z4 -H- — z). Furthermore, since Z4 becomes unbounded as ?/ — ?■ ?/o, eq. ( 1T8|) 
implies that the submatrix — Z4 is the conditional impedance with the auxiliary condition 
of rigid (infinite) impedance at y = yo, an observation that is verified by eqs. ([9]) 2 and 
dUi. 



2 Cylindrically anisotropic elastic solids 
2.1 Equations in cylindrical coordinates 

The dynamic equilibrium equations for a linearly elastic material when expressed in cylin- 
drical coordinates are jl^ 

/O -1 0\ 

r-^(rt^),^ + r-^(te,0 + Kte) = pii with K = 1 . (19) 

\0 0/ 

Here p = p(x) is the mass density, u = u(x, t) the displacement, and the traction vectors 
tj = tj(x, t), i = r,6,z, are defined by the orthonormal basis vectors {6^,69, e^} of 
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the cylindrical coordinates {r,6,z} according to tj = e^a {i = r,6,z), where <t(x, t) is 
the stress, and a comma denotes partial differentiation. With the same basis vectors, 
and assuming the summation convention on repeated indices, the elements of stress are 
aij = CijkiSki where e = |(Vu + Vu-^) is the strain, Cijki = Cijkii'x.) are elements of the 
fourth order (anisotropic) elastic stiffness tensor, and T denotes transpose. The traction 
vectors are 
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where, in the notation of the matrix (ab) has components {ah)-f^ = aiCijkibi for arbitrary 
vectors a and b. The explicit form of the various matrices is apparent with the use of 
Voigt's notation Cijki — Cap (a, /? G {1,2,..., 6}) 
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C16 C12 Ci4 
C66 C26 C46 
vC56 C25 C45 



(20) 



2.2 Cylindrically anisotropic materials 

The concept of cylindrical anisotropy, which appare ntly originated with Jean Claude 
Saint- Venant, and has been elaborated by Lekhnitskii [2^, demands the angular indepen- 
dence of material constants in the cylindrical coordinates, but admits their dependence on 
r and as well on z. We consider materials with no axial dependence whose density and the 
elasticity tensor may depend upon r, i.e. p = p{r) and Cijki = Cijkiir) \/i,j,k,l G r,9,z. 
We seek solutions in the form of time-harmonic cylindrical waves as 

where ?T, = 0,l,2,...is the circumferential number. 

The dependence of the displacement and traction on the single spatial coordinate r 
allows the elastodynamic equations to be reduced to the canonical form of eq. ([T]) |l2l ]: 



^^W(r) = ^G(r)i7(")(r), (22) 



where rj^''^^ is a 6 x 1 vector 



V^'^Hr) = (^Sl5|"|^) , with V(")(r) = ^rTW(r), (23) 
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and the 6x6 system matrix G is defined by 



iG{r) = go(r) + rgi(r) + r^g2{r) 
The individual 6x6 matrices are 





So — I ^"isl I ' gl — "^^2 I ^^3} inr 1 ; §2 — ( . {3} r> 

' ^To -go / V^gi gi / V*g2 

with the 3x3 matrices 

g^'^ = gi'^ + ^Krgl'\ g^'^ = gl'\ g^'^ = g^'^ + ^hrgi'^+r'gi'\ 
The constituent 3x3 matrices are 

gj^> = -Q-iR, g« = -Q-^P, 

si'' = -Q ' = gj^^^, gP^ = P^Q"^R - S - (P-Q-^R - S)+ = -gf >^ 

gj^> = T - R+Q-iR = gf >^ gf = A:,^(M - P^Q-^P) - pujH = gf 

where 

R = Rk, S = «S, T = /c+T/c = T+, K = K + mI = -/c+. 

The matrices go^-*^ and g^^^ are negative definite and positive semi-definite, respectively, for 
real-valued and positive definite elastic moduli. Note that the nth order modal solution 
rj^'^\r) is a function of the radial coordinate, but it is also an implicit function of the 
frequency u and the axial wavenumber k^, which dependence is here kept tacit. In the 
same manner, the dependence of G(r) upon n, u and is understood. 

The superscript '•"^ is omitted henceforth, with the exception of the specific cases 
n = 0, 1, as required. 



2.3 Cylindrical elasticity in the general context 

The results of the previous subsection, particularly eqs. (!22l) and (!23l) . show that the 
cylindrically anisotropic system of azimuthal order n is a special case of the formulation 
of Section [1] generally with m = 3, and {y, U, V, Q} — )■ {r, U , irT, ir~^G}. The physical 
restrictions required for the Hermiticity condition are real-valued y, fc^ and material 
constants (more precisely, Hermitian elastic moduli c^/? = c^^ suffice Under these 

conditions the 6x6 matrix G(r) displays the symmetry 

G = TG+T. (24) 

The 6x6 matricant M(r, ro) is the solution of the initial value problem 
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The condition that r and tq are strictly positive is important since the case of zero radial 
coordinate needs to be handled separately, which is discussed at length below. Note that 
we do not specify whether r or tq is the greater or lesser of the two radii. The matricant 
allows us to express the state vector rj (r) of partial modes in a cylinder as 

r]{r) = M (r, Tq) T]{ro), r, Tq ^ 0. (26) 

The pointwise elastodynamic energy balance is dS/dt + divV = where S is the 
energy density per unit volume and V the energy flux vector. The pertinent form of ([2]) 2 
for cylindrical elasticity is divP = r~^d{rPr)/dr = where Pr = {V)t ■ e,,. is the time 
averaged radial component for azimuthal mode n, 

P,(r) = -£-f/+(r)Tr/(r), (27) 

which together with the system equation fl22|) implies the symmetry fl24|) for G (see also 
0)- 

The conditional impedance matrix z relates traction and displacement at a particular 
value of r, but specifically r 7^ 0, according to eq. IQ. The point r = requires a separate 
discussion, and indeed a newly defined impedance, introduced in the next section. For the 
moment we note that z(r) is contingent upon the definition of the (one-point) impedance 
at some radial coordinate, say z(ro) = Zq. The traction at other values of r is then 
unambiguously related to the local displacement by either the matricant or the two-point 
impedance matrices, using equation or (fT7|) . By rewriting eq. fl27p we see that the 
conditional impedance determines the pointwise flux, 

p(r) = -— Im{U+(r)z(r)U(r)|, (28) 
2r ^ ^ 

which is zero for all U(r) only if z is Hermitian. This in turn is the case only if z(ro) = zq 
is Hermitian, i.e., if there is no flux across the surface r = tq. On the other hand, the 
6x6 two-point impedance matrix Z(r2,ri) of eq. (fT2l) defines the global energy flow into 
or out of the finite region between the two radial coordinates ri < r2. Let E{t) be the 
total energy in the shell cross-section per unit length of the cylinder for azimuthal mode 
n. Its increment over one period of time harmonic motion is 

-2,f - -2.^.n,{ (u(-))^z(,..,..) (j;<->) }. (29, 

which is identically zero for real u, kz and Hermitian parameters, i.e., when Z is Hermitian. 
If the material in the slab is lossy then Im(Z — Z"*") should be positive definite in order 
that E is not increasing with time. 

The differential Riccati equation satisfied by the one-point impedance matrix follows 
from ( 1T0|) as 

+ zg{i} + g{i}+z + zg{2}2 + g{3} ^ 0. (30) 

dr 

The initial value problem for z(r) is therefore 
r— - [z+ (R + zA;,rP)^]Q-^[z + R + zA;,rP] +B(r) = 0, r > 0, z(ro) = zq, (31) 
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where 

B(r) = f + ik,r{S+ - S) + r^k^M - puH) = B+(r). (32) 

Equation fl3T]) shows the exphcit dependence upon u, and the elastic moduh. The 
exclusion of the distinguished point r = at the cylinder centre is addressed next. 




Figure 1: Three types of cylindrical structures defined hj ri < r < r^. the annulus 
(0 < ri < r2 < cxo), the solid (ji = 0), and the exterior region (r2 = oo). 



3 Wave impedance matrices for cylinders 

In this section we describe typical uses of impedance matrices, and in the process introduce 
the solid-cylinder impedance Z(r) and the radiation impedance Zi.ad(^)- We consider the 
three distinct configurations depicted schematically in Figure [TJ 

3.1 Solid-cylinder impedance matrix Z(r) 

A solid cylinder, by definition, is one that includes the axis r = 0. A new impedance matrix 
is introduced to handle this situation. The solid-cylinder impedance Z(r) is defined in 
the usual manner by its property of relating the traction and displacement 3— vectors of 
f l2T|) according to 

V(r) = -iZ(r)U(r), r > 0, (33) 

although this is not a conditional impedance matrix because of the absence of an auxiliary 
impedance condition at some other coordinate. Instead, the solidity of the cylinder at 
r = dictates the character of Z(r) (and one could argue that it is 'conditional' in that 
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sense). The limiting value of the solid impedance at r = plays a crucial role, and we 
accordingly define the central-impedance matrix: 

Zo = Z(0). (34) 

The properties of the central impedance are discussed in detail in ^ after we develop 
methods for finding the solid-cylinder impedance matrix in ^ 

As an example application consider the task of finding the dispersion equation for 
guided waves of frequency u and wavenumber k^. We suppose, quite generally, an interface 
condition on the level surface r = r2 of the form 

V(r2) = -?Z2U(r2), (35) 

where Z2 is considered as given. It could be zero (traction free condition), infinite (rigid 
boundary), or it could be defined by some surrounding material, whether finite or infinite 
in extent. For instance, if the solid cylinder is surrounded by a shell of cylindrically 
anisotropic material in lubricated contact at r = r2 and free at r = rs > r2, then Z2 = 
Zii{r2)erej where z(r2) is the conditional impedance with the auxiliary condition z(r3) = 
0. Assuming fl35|l describes the condition at the outer surface, the desired dispersion 
equation is 

det (Z(r2) - Z2) = 0. (36) 

It is instructive to compare ( 136|) with the dispersion equation for a (possibly functionally 
graded) layer y e [0,1/2] with the traction- free surface ?/ = on a homogeneous substrate 
y > 1/2 which may be written in the form jUj 

det (z (2/2) - Z2) = 0, 

where Z2 is a (constant) impedance of the substrate and z (7/2) is the conditional impedance 
of the layer satisfying the reference condition z (0) = 0. If the surrounding material beyond 
a rigid (say) interface r = r2 is infinite, then there are Stoneley-like waves defined by the 
dispersion equation f p5|l with Z2 = Zrad (^^2) , where Zrad {r) is the radiation impedance 
discussed in ^ 

The solid-cylinder impedance also provides a means to compute the modal displace- 
ment vector U(r) for all < r < r2 if the dispersion equation f l36|) is satisfied. By analogy 
with the case of an annulus j2Tj [. the unnormalized displacement follows from the system 
equation fl22|) and the definition of Z(r) in f l33|) as the solution of the initial value problem 

r^ + Q"'(R + ^A;,P + Z)U(r) = 0, 0<r <r2; U(r2) = Uq, (37) 
dr ^ 

where Uq is the null vector of the surface impedance condition, (Z(r2) — Z2)Uo = 0. Note 
that the solution of fl37|l remains well behaved even as the matricant solution U(r) = 
(Mi(r,r2) -iM2(r,r2)Z(r2))U( r2) is numerically unstable, see §4.2.21 

3.2 Impedance matrices for cylinders of infinite radius 

Consider a cylinder extending to infinity in the radial direction, with inner surface at 
r = ri, see Figure [1] A wave incident from r > ri results in a total field that can be 
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expanded in terms of partial waves of the form fl^T]) . The amphtude of the n*^ azimuthal 
mode is 

U(r) = Ui,,(r) + U,,at(r), r > n, (38) 

where the scattered amphtude Ugcatl'") satisfies a radiation condition at r — > oo. This in 
turn requires that the following condition prevails on the interface: 

Vscat(r) = -iZrad(r)Uscat(?^), T = n, (39) 

where the radiation impedance matrix Zi.ad('") is defined by the radiation conditions, see 
^ The scattered field is then uniquely determined by the condition at r = ri, which we 
assume is of the generalized form (135|) with prescribed interface impedance zi. Then, 

- Zinc(ri)Uinc(ri) - Zrad(?'l)Uscat('"l) = -Zl(Uinc(ri) + Uscat(n)), 

where Zinc(?") is the impedance of the incident wave, which follows directly from the 
equations of motion. The scattered amplitude on the interface is therefore 

Uscat(ri) = (zi - Zrad(ri))"^(zi - Zinc(ri)) Uinc(ri), (40) 

which provides the initial condition to determine the entire scattered field in r > ri. Fur- 
ther details on the radiation impedance matrix are provided in §6l including its asymptotic 
properties for large r. 

3.3 An annulus of finite thickness 

The case of the annulus < ri < r < r2 fits readily into the general theory. Again 
consider the task of finding the dispersion equation for guided waves, which may be found 
by simultaneous satisfaction of the conditions on the two radial surfaces. Suppose the 
conditions are both of the generalized form V(rj) = — zzjU(rj) (j = 1,2) where zj [j = 
1,2) are known quantities. The conditional impedance z(r) is determined (numerically) 
by integrating (1251) from (say) r = ri with initial condition z(ri) = zi to give 

z(r) = i(M3(r,ri) -iM4(r,ri)zi)(Mi(r,ri) -iM2(r,ri)zi)"\ (41) 
The interface condition at r = r2 requires that 

- Z2U(r2) = -z(r2)U(r2), 
which implies the dispersion equation 

det{z(M3(r2,ri) -zM4(r2,ri)zi)(Mi(r2,ri) -zM2(r2,ri)zi)~' -Z2} = 0. (42) 

Variants on this equation may be obtained using the two-point impedance instead of the 
matricant. Thus, from eq. f|T6l) we have the equivalent condition 

det {z2 + Z4(r2, ri) + Z3(r2, n) (zi - Zi(r2, n)) "'Z2(r2, n) } = 0. 

The examples considered in this section illustrate the usefulness of wave impedance 
matrices for cylinders of finite and infinite radial extent. Solutions to problems of practical 
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concern can be formulated concisely in terms of impedance matrices, such as the dispersion 
equation for guided waves, and the scattering of waves from a cylindrical region. Cal- 
culation of the impedance matrices is relatively straightforward using the matricant or 
two-point impedance matrices (see [lH), but only as long as the points r = or r = oo are 
not involved; otherwise the solid cylinder impedance and/or radiation impedance matri- 
ces are required. Determination of the solid-cylinder impedance matrix Z(r) is discussed 
next. 



4 The solid impedance matrix 

In this section we develop methods to calculate the solid-cylinder impedance matrix for 
a radially inhomogeneous cylindrically anisotropic cylinder with material at r = 0. Two 
principle approaches are considered: a semi-explicit solution as a Frobenius series, and an 
implicit solution in terms of a differential Riccati equation. 

Unlike the conditional impedance which can be determined directly from the matricant 
M along with the prescribed reference value, the matricant is not of direct use here because 
of its divergence at r = 0. This introduces the need to identify 'physical' and 'nonphysical' 
constituents of the solution near r = 0, which is performed explicitly for the Frobenius 
solution. In the Riccati approach the displacement and traction fields are not considered 
explicitly and the divergence at r = is taken care of by the initial value of the impedance. 
The Frobenius solution is considered first. 



4.1 Frobenius expansion 

We take advantage of the fact that the fundamental solution can formally be written 
in terms of a Frobenius series, which is an explicit one-point solution valid at any r 
(including r = 0). As a result, the Frobenius-series approach provides a constructive 
definition of Z(r). The Frobenius series solution can be obtained via a recursive procedure 
with the number of numerically required terms increasing with r. Before we present the 
formal solution for Z(r) we review and develop some properties of the Frobenius series 
for cylindrically anisotropic materials, following the analysis in [20| . 

4.1.1 Background material 

The Frobenius solution is based on the integral matrix solution M {y) = \\tji, ...,t}q\\ of 
Eq. (1221) . which can always be defined through the Frobenius series for any r > 0. The 
pivotal role in constructing this series belongs to the eigenspectrum of the 6x6 matrix 
go (0) with the symmetry 

go = -Tg+T, (43) 

which follows from flMl) . Denote the eigenvalues and eigenvectors of go (0) by Aq, and 

7oa = {sLa,ia) (o = 1,...,6), and introduce the matrix To = UToi, •••^Toell- Barring 
extraordinary exceptions, if n > 1 then (i) no two eigenvalues Aq, of go (0) differ by an 
integer, and (ii) all Aq, normally are distinct (and nonzero). Let us first consider this case 
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n > 1, otherwise see §4.1.2[ By virtue of (i), the integral matrix may be formulated as 

oo 

j\A (r) = D (r) ToC, D (r) = I + D^r"^, (44) 

m=l 

where C is the Jordan form of the matrix r^"^^^ which is diagonal when (ii) holds, and 
D (r) is defined recursively through G (r) |20|, Eqs. (9)-(13)]. 

The arguments underlying eq. ([2]) imply that the matrix jV+ (r) TjV (r) is a constant 
independent of r, and according to eq. (1271) this matrix defines the flux properties of the 
constituents r/i, . . . , r/e, of jV (r). For the present purposes we wish to split them into a 
pair of triplets: a physical set (a = 1,2,3) and a nonphysical triplet (a = 4,5,6), where 
the only non-zero flux interactions occur between a and a + 3 (a = 1, 2, 3), thus ensuring 
the crucial property that jV^ (r) TjV (r) has nonzero elements confined to the main di- 
agonal of the off-diagonal blocks. The partitioning is accomplished through appropriate 
arrangement of the eigenspectrum of go (0) as 

K = -K+3^ ReA, >0, a = 1,2, 3, (45) 

|20l . Eq. (44)]. Combining Eqs. ( l43l) and ( 145|) and adopting the normalization 7QQ,T7oa+3 = 
1 yields the orthogonality/completeness relation for the eigenvectors in the form 

r+TTo = T. (46) 

It follows from Eqs. through that jV+ (0) TjV (0) = T and hence the flux matrix 
at r is T, 

(r) TAA (r) = T ( ^ jV (r) TM^ (r) = T) . (47) 

Note that yields D+TD = T. 

In order to further clarify the structure of J\f we represent the 6x6 matrices D, Tq 
and C in terms of 3 x 3 submatrices, 

/Di D2\ /Ai AaX /diag(r^'^) \ 

D(r)= hro= 'C= , (48) 

\D3 D4/ V^i L2/ \ diag(r-^^)/ 

where a = 1, 2, 3 and C is diagonal for n > 1. The integral matrix consequently has block 
structure 

/Ui UaX /Di DaX /Aidiag (r^«) Aadiag (r-^S)\ 
Af{r)= ] = [ . (49) 

\Vi V2/ \D3 D4/ \Lidiag(r^") Lsdiag (r"^-) / 

Note in particular that the integral matrix A/" (r) consists of two distinct 6x3 matrices, 

Vi,V2,V3\\ = D ||7oi,7o2,7o3|| diag (r^") , 

(50) 

V4,V5,V6\\ = D ||7o4,7o5,7o6|| diag (r"-^-) , 
14 





the former with the columns f/a (r) tending to zero at r — )■ 0, and the latter with columns 
Va+3 diverging at r — )■ 0. The block structure of eqs. (jlH]) and fll7|) is 





(51) 



jV+ (r)Tj\f{ 



'^t V+\ /Vi 

SJt ytJ \ui 




The latter explicitly shows that the normal energy flux of the displacement-traction wave 
field rj (r) comprising an arbitrary superposition of either the three modes rj^ (r) or three 
modes t/q+s (r) with a = 1, 2, 3 is zero at any r. This specific arrangement of Af may be 
interpreted as the generalization of the isotropic case with solutions cast in terms of the 
cylinder functions J„ and —iYn, corresponding to the physical and nonphysical triplets 
respectively, each of which yields zero flux individually. This partitioning will be crucial 
in developing an explicit solution for the solid impedance matrix. 



4.1.2 Overview of the cases n = and n = 1 

Let us return to the two assumptions made above which are that (i) no two eigenvalues 
Aq of go (0) differ by an integer and (ii) all Aq are distinct, hence go (0) is semisimple 
(diagonalizable). Violating (i) invalidates the relatively simple form (jH]) of the Frobe- 
nius fundamental solution to the governing equation, see |23j. Violation of (ii), or more 
precisely, the occurrence of degenerate Aq that makes go (0) non-semisimple, alters the 
orthogonality/completeness relations and the composition of J\f given above for tt, > 1. 
The cases affected are n = (axisymmetric modes) and n = 1 (lowest-order fiexural 
modes): specifically, the property (i) does not hold for n = 0, and the property (ii) does 
not hold for both n = and n = 1. From a physical point of view, the cases n = 0, 1 
stand out because they are related to the rigid-body motions producing zero stresses |20l . 
Eq. (19)]. Note also that go (0) admits a zero eigenvalue iff n = 0, 1 [20|, Eq. (80)3] and 

that A*-"'^-* = is always a double eigenvalue rendering go'^'^^ (0) non-semisimple. 

Consider the axisymmetric case n = 0. The six eigenvalues Xa'^ of gg'^'' (0) are Xa^ = 
{0, 0, ±1, ±/t} , where k = 1 for trigonal or tetragonal symmetry with Ciq = [13, Eqs. 

(3.12), (3.13)]. It is seen that, whatever the symmetry, the set of Xa^ includes pairs 
different by an integer. As a result, the integral matrix J\f (r) is now defined through 

gl^^ (0) in a rather intricate form elucidated in jlol, Eqs. (A2), (A. 4)]. This observation 
is essential for treating inhomogeneous and low-symmetry homogeneous cylinders. At 
the same time, if the cylinder is homogeneous and has orthorhombic or higher symmetry 
with the exception of trigonal and tetragonal with cig = 0, then J\f (r) decouples into the 
solutions described by Bessel functions and/or by a simple Frobenius form ( H^ll PI. 

"'^ Orthorhombic or higher symmetry enables uncouphng of the pair of torsional modes described by the 
Bessel solutions stemming from A'°' = ±1. The four sagittal modes are associated with A'-^-' = {0, 0, ±k} , 
where k ^ I for symmetry lower than the trigonal or tetragonal with cie = 0. When k = 1, so that the 
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Consider the case n = 1. The matrix (0) has a doubly degenerate eigenvalue 

A*-^^ = which makes non-semisimple j20|, Eq. (36)]. This does not preclude taking 
jV(r) in the form (jH]) but the matrix C is now not diagonal. As a result, the triplet 

q; = 1, 2, 3 of physical modes (with one of the modes r]a^ associated with A*-^^ = 0) retains 
its form ( l50|) i . whereas the nonphysical triplet a = 4, 5,6 is no longer of the form (!50|) '> 
due to one of its modes involving both ei gen vectors, the proper and the generalized ones 
7^-^) and 7*^^\ associated with A*^^) = ]20l . Eqs. (51), (61)]. It is thus evident that 
the physical modes satisfy the same orthogonality/completeness relations as for n > 1; 
moreover, subject to the optional condition ^(^)+Tf^^^ = 0, the nonphysical modes may 
be shown to do so as well. The relations ( 146|) and ( 1471) for the case n = 1 are accordingly 
modified into a slightly different form 

r+TFo = E, J\r+ (r) TAT (r) = E (n = 1) , (52) 

which differs from ( H6l) and ( H71) only in the replacement of the right-hand matrix T by E, 
whose nonzero elements are also confined to the main diagonal of the off-diagonal blocks 
but they cannot now be all normalized to 1 jlol, Eq. (49)]. 

The overall conclusion is that both cases n = and n = 1 preserve the partition- 
ing of the six linear independent Frobenius solutions rfa (r) = (Uq,, Vq,)^ within M {r) 
(a = 1,...,6) into the physical and nonphysical triplets a = 1,2,3 and a = 4,5,6. The 
partitioning is based on fjlSj) supplemented by including the (double) eigenvalue A^^^'^^ = 0. 
The vectors Uq, (r) and Vq, (r) are certainly regular at r — )■ for both n = 0,1, although 
the limiting trend for n = is not of the form that results from (jH]), see (20|, Eq. (A4)]. 
Equations f l50|) and f l5T|) . which are valid for any n > 0, enable treating the solid cylinder 
impedance Z (r) for n = 1 on the same grounds as for the 'ordinary' case n > 1. The 
impedance Z (r) for n = needs special attention because the case n = may not satisfy 
dH]) . We are now ready to derive the explicit form of the solid-cylinder impedance for all 
n. 

4.2 Explicit solution of the solid-cylinder impedance 
4.2.1 The solid-cylinder impedance for arbitrary n 

The definition ( 1331) of the solid- cylinder impedance Z (r) tacitly assumes U(r) and V(r) 
are regular function of r. This is always so for rf (r) = (U, V)""" comprising an arbitrary 
superposition of, specifically, the physical Frobenius modes r)a (r) = (Uq, V^) which 
satisfy eq. f H5|) supplemented by the option A*-"'^-* = for = 0, 1 (see §4.1.2p . Thus the 
solid-cylinder impedance may be defined by any of the equivalent expressions 

Va{r) = -iZ{r)Va{r) (a = 1,2, 3) ^ V (r) = -zZ (r) U (r) ^ Z (r) = zVi (r) U^^ (r) . 

(53) 

This yields a finite value if det Ui (r) 7^ 0, otherwise the impedance is associated with a 
'rigid' condition at r (conversely, the determinant of its inverse - the admittance matrix 

above quartet of A*-*^^ involves pairs with an integer difference, tlic sagittal problem admits explicit Bessel 
solutions for the isotropic or transverse isotropic symmetry due to uncoupling of potentials. Note that 
double eigenvalues A*-"^ = ±1 at k = 1 do not bring non-diagonal blocks into the Jordan form of g'^^ (0) . 
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- is zero). The occurrence of infinities is in no way anomalous but rather a natural 
consequence of the definition of the impedance matrix. 

Consider first n > 0. Based on the definition flS^ and the representation flSD]) i for the 

3x3 matrices Ui and Vi, we obtain an alternative form for the solid-cylinder impedance, 
Z(r) = 2(D3(r)-zD4(r)Zo)(Di(r)-iD2(r)Zo)"' where Zq = zLiA^^ (54) 

Hermiticity of the solid-cylinder impedance follows from f lSTjl ? and fl52l) 9 which imply that 
U|Vi + V^Ui = -ilJi (Z - Z+) Ui = 0, whence 

Z (r) = Z+ (r) . (55) 

The expression (l54l) is reminiscent of the representation of the conditional impedance, 
e.g. eq. (SI]), except that the role of the two-point matricant M(r, tq) is replaced by 
D(r). Note that det Ai ^ may be deduced from the integral representation of Zq, see 
eq. ( 179|) i ■ by reasoning similar to jH: that if two of the eigenvectors are parallel, say 
ai and sl2, then so are the traction vectors Iq, = —iZoSLa counter to the assumed linear 
independence of 71 and 72. 

Now consider n = 0. Violation of eq. (l44l) for A^^^-* (r) invalidates the definition ( l54|) '? 

for the central impedance Zq'^''. At the same time, Z^^ can readily be found by means of a 
direct derivation given in §5.4.21 specifically eq. ( l90l) . which is clearly Hermitian regardless 
of anisotropy. Consequently Z^'^^ (r) is Hermitian for any r due to the self-adjoint property 
of the differential Riccati equation of which Z^^^ (r) is the unique physical solution (see 



4.2.2 The link between the solid-cylinder and the conditional impedances 

It should be evident from the previous discussion that Z (r) can formally be defined as 
the conditional impedance z (r) with initial value at tq — > 0. Assume for brevity that 
n > 1, then using the representation M(y,?/o) = {y)J^~^ iVo) for the matricant and eqs. 
(iZD, (USD, we have 



M(rr) -4 {^'^''^ U2 (r) \ /V2+ (ro) U+(ro)A Ui(r)V2+(ro) Ui(r)U+(ro) 

^ ' °^o-.o I Vi(r) V2(r)i V ) lVi(r)V2+(ro) Vi(r)U+(ro)^ ^ 

This illustrates that even though the matricant M (r, tq) diverges at tq — > 0, as expected, 
it provides the correct limit 

z(r) = z(M3-«M4z(ro))(Mi-zM2z(ro))~^ ^ zVi (r) U^^ (r) = Z (r) . (56) 

ro^-O 

Formal consistency requires the limiting value of z (ro) at ro — be set equal to Zq. 
However, the definition (156|) of Z (r) is actually of no value for practical calculations 
because of the divergence of M (r, ro) at ro — 0. 

At the same time, in the limit as r — > the conditional impedance z (r) with any 
initial value z (ro) such that |z(ro) — Z(ro)| > 0, should tend to the nonphysical central 
impedance Z^p (0) = iL2A2 ^. Similarly to ( 156|) . 

M(rr) ^ ^"^A ( ^'''^ ^'''^ ] = l^' ^"^^ ^"^'^ ^""^ ^"""^ 

^ ' °^-.o I V2 (r) i I V+ (ro) U+ (ro) / I V2 (r) V+ (ro) V2 (r) U+ (ro) 
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Hence, from fl4ip . 



(r) = z (M3 - zM4Z (ro)) (Mi - zMsZ (ro))' 



r->0 



V2 (r) U2 ^ (r) 



r-i-O 



zL2A, 



-1 

2 ■ 



(57) 

If z(ro) is precisely the solid impedance at ro then f l57|) i reproduces the solid impedance, 
z (r) = Z (r) for r > 0. But the limit at r = 0, formally Z (0) = Zq, cannot be achieved 
in practice, a reflection of the fact that the matricant based solution ([7]) in cylindrical 
coordinates is uniquely ill-posed at this point (see also §4.4p . 



4.3 Riccati equation solution 

An alternative to the Frobenius approach is to consider Z(r) as a solution of the differential 
Riccati equation with initial value Zq extended to the case when the initial value occurs 
at r = 0. The solid- cylinder impedance is then the solution of the initial value problem, 

= [Z+ {R + ik,rPy]Q-^[Z + R + ik,rP] -B{r), r > 0; Z(0) = Zq, (58) 

where B(r) is defined in eq. (132|) . The central-impedance matrix Zq, as discussed in 
the previous subsection, is defined by the eigenvectors of go(0), see (154|) 9. Alternatively, 
noting that a nonphysical singularity is introduced unless the right hand side of (!58ll i 
vanishes at r = 0, we conclude that the central impedance must satisfy the algebraic 
Riccati equation 

(Zo + R+) Qo ^ (Zo + Ro) - To = 0. (59) 

While it is expected that the solution Z(r) is well behaved in some finite neighbourhood 
of r = 0, the Riccati solution will inevitably develop singularities. These are associated 
with guided waves of a cylinder of radius r with clamped surface (zero displacement 

condition). For given u and k^, the singularities occur at values of r such that det Ui (r) = 
(see eq. (l53l)). Thus, one can integrate the differential Riccati equation only as far as 
the first singularity at (say) r = r*. The problem is evident from the example of the out- 
of-plane impedance derived in f lllOp -). Zz{r,0) = —CMk2rJ^{k2r)/Jn{k2r), which blows 
up when k2r is a zero of the Bessel function J„. The effect of singularities may be 
circumvented in practice by integrating the Riccati equation to some finite r short of the 
first singularity and then to switch to some other solution method that is regular in the 
vicinity of r = r*. One approach jsj is to consider the admittance (inverse of impedance) 
Y(r) = Z~^(r) which will be well behaved at r = r^,. Its differential Riccati equation, 
which is easily found from ( l58l) . can therefore be integrated without incident through 
the singularity at r = r*, but the admittance then has its own singularities at positions 
different from those of the impedance, so in general this approach requires switching 
back and forth between two Riccati equations. While certainly feasible, the procedure is 
complicated by the fact that one does not know the singularities a priori. Note that the 
admittance Riccati equation is not suitable for starting at r = because as discussed in 
the next section det Zq = and hence Yq = Y(0) is undefined for n = 0, 1. 

A more practical approach to deal with the unavoidable singularity problem is to use 
the Riccati solution to generate initial conditions for the full 6x6 system at r = ri < r*. 
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with which one can integrate (again numerically) to arbitrary r > ri using fHTl) . In 
practice one only needs to solve for a 6 x 3 matrix fl{r), satisfying 

A^2(r) = lG(r)f2(r), r > n; l^(rO = (_^2(ri)) ' 

Although Q{r) does not describe the complete wave field it is sufficient to determine the 
impedance, since 

= {n'^ = (-zzVi)) ^ ^^"^ = ^n,{r)n^\r), (61) 

for r > ri. The value of r at which one switches from the differential Riccati equation to 
the matricant based solution is a free parameter, and arbitrary as long as it is lies below 
the first singularity in the impedance. This can be estimated from the separable solutions 

1 /2 

in ^as r^(uj'^s'^^^ — k^) ~ 1 where Smax is the largest plane wave slowness at r = 0. 
4.4 Discussion 

We have described two principal ways for finding the solid-cylinder impedance Z (r). The 
Frobenius series method is summarized in Eqs. f l53|) and fl5^ . Taken together these 
equations provide a basis for calculating the solid-cylinder impedance for n > and 
arbitrary r via the Frobenius series solution. The Riccati equation method determines 
Z (r) for arbitrary n by integrating the differential Riccati equation fl58|) subject to an 
initial condition defined by the central impedance Zq. The Riccati approach is strictly 
valid only for r less than the first singularity of the solid-cylinder impedance. 

The initial value Zq can be evaluated from ( 15^ 9 or by other methods discussed in ^ 

For n = the form of Zq''^ is explicit (eq. ( I90l) below) and Z^^^ (r) may be determined by, 
for instance, integration of the Riccati equation discussed in §4.31 The physical solution 
to the initial value Riccati equation can be continued through and beyond the first and 
subsequent singularities by using the matricant solution to generate Z(r) as a conditional 
impedance. Strictly speaking the practical value of the Riccati method is confined to 
the neighborhood of r = 0. The differential Riccati equation provides a regularization 
of the system of equations ( l22l) which are singular at r = 0. Once this singularity has 
been taken care of, there is no need to use the Riccati equation, particularly since the 
Riccati equation has its own singularities - in fact an infinite number of them. Note that 
satisfaction of the algebraic Riccati equation (!59l) is essential to ensure regularization of 
the initial value problem ( l58|) at r = 0. The differential Riccati equation cannot generally 
recover the central impedance Zq by 'backward' integration to r = from some initial 
To > 0, because the system possesses the same ill-posed property observed with respect 
to eq. f l57|l . in this case associated with the fact that the nonphysical central impedance 
Znp (0) (= A^^ for n > 1) also solves (159|) . 

Both the Frobenius and Riccati methods generate an Hermitian solid cylinder impedance. 
Hermiticity of Z(r) is a consequence of the fact that it is built from the triplet of physical 
modes which produce zero normal fluxes both of their own and due to their cross-coupling. 
Note that the nonphysical impedance Znp (r) = zV2U2^^ is Hermitian as well, which is 
similar to the case of a half-space; however, the physical and nonphysical impedances 
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of a cylinder are generally no longer negative transpose of each other as they are for a 
half-space. For n > 1, the two impedances are related by 

Z(r)-Z,p (r) = z(UiU+)-\ 

with normalized Ui, U2 on the right hand side, as follows from f H7|l 9 and fl^Tlj o. The 
Hermitian nature of Z(r), r > can also be viewed as a consequence of the fact that it 
solves the Riccati equation ( !58l) with an Hermitian initial value, Zq. It is also noteworthy 
that neither the definition ( !53l) of Z (r) nor its Hermitian property requires any specific 
normalization of the eigenvectors 7^ of go (0) once they have been ordered into physical 
and nonphysical triplets. 

While the solid-cylinder impedance is quite distinct in nature, it is in a certain sense, 
a conditional 'one-point' impedance, for it depends on the initial condition at r = 0. 
However, there is another, more essential, aspect, which actually sets Z(r) apart from 
the two-point impedance Z (r, tq) and the general conditional impedance z(r). It is that 
Z (r, tq) and z(r) involve all 6 linear independent partial solutions, whereas Z(r) involves 
only half of them and discards the other half on the basis of certain partitioning at 
r = (physical/nonphysical). As a result, the Hermiticity of Z (r, tq) and z(r) and the 
Hermiticity of Z(r) have different origins. Hermiticity of both Z (r, tq) and z(r) (the latter 
subject to Hermiticity of the initial condition) follows from divP = while Hermiticity 
of Z (r) is a consequence of = 0. 



5 Properties of the central impedance matrix Zq 

The central-impedance matrix depends only on the elastic moduli (and n) and is simpler 
than the solid- cylinder impedance, its continuation away from r = 0. At the same time, 
the value of the central-impedance is required a priori in order to calculate Z(r) using 
the Riccati equation (!58|) . In this section we describe some properties of Zq, develop 
procedures for its determination, and consider its behavior for large n. 



5.1 Integral formula for Zq with n > 1 

5.1.1 The Lothe-Barnett integral method using the matrix sign function 

The surface impedance matrix Z(w), v =velocity, plays a central part in the theory of 
surface waves in an elastic ho mog eneous half-space. It was first identified in that context 



by Ingebrigtsen and Tonning 2^ and subsequently developed as a crucial ingredient for 
pro ving the uniqueness and the existence conditions for surface waves [H, Q], see also 
[281 . 129I]. The central impedance matrix Zq of a cylinder has a close relationship to the 
static [v = 0) surface impedance matrix. The similarity allows us to use some of the 
considerable array of results for the latter. Here we draw directly on the integral formalism 
for the surface impedance of a half-space first outlined by Barnett and Lothe |30l. l3ll ] and 
then presented in full by Lothe and Barnett [l| and Chadwick and Smith |28l |. In this 
subsection we show how this formalism can be modified to describe the cylinder central 
impedance Zq for n > 1, and we discuss the exceptional cases n = 0,1. 

Assume n > 1, so that the eigenspectrum of go(0) lies on either side of the imaginary 
axis in accordance with eq. (145|1 . The restriction to n > 1 will be clarified below. The 
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matrix sign function is then uniquely defined 

signgo(0) = go(0)(go2(0))^'/', (62) 

where the principle branch of the square root function with branch cut on the negative 
real axis is understood; z = [z^Y^"^ sign 2; with sign 2; = +1(— 1) if Ylez > 0(< 0). As a 
result the sign matrix satisfies 

( signgo(0))7a = ±7a for ReA„ ^ 0. (63) 

The matrix sign function was first introduced by Roberts [S^j as a means of solving 
algebraic Riccati equations, and has become a standard matrix function, see jsH, 113] for 
reviews; the simple relation fl62|) was first noted in jlH. Using the spectral decomposition 
defined by the matrix of eigenvectors Fg yields 

go(0) = roAoFo ' signgo(O) = ro(signAo)roi with 

^diag(A,) \ ^\ 

Ao= I \ ^ signAo= . (64) 

diag(-A;)/ \0 -I J 



The explicit structure of the sign matrix follows from the normalization condition (Il6 
and the sub matrices defined in ( HHII q. 



sign go(0) = To ( sign Aq) TF+T = (^^ with 
S = 2A1L+ -1 = 1- 2A2L+, H = -2iAiA+ = H+, B = -2iLiL+ = B+. (65) 
Additional relations are obtained from the involutory property of the sign matrix function, 

(signgo(O))^ = I6 S^-HB = I, SH=(SH)+, BS = (BS)+. (66) 

The connection with Barnett and Lothe's theory is established via the integral expres- 
sion for the matrix sign function 32|, |34 



signgo(O) = -go(0) / dt {t'l + g^(0))"\ (67) 
^ Jo 



A simple change of integration variable and separation into partial fractions yields sign go(0) 
as an averaged matrix. 



7r 

signgo(0) = i I d0gf = (gf >, (6J 







where 



Q — (cos0I — 2sin0go(O)) ^ ( COS0 go(0) — i sin0 1). (69) 



— I L,un 1/7 J. — t niii 1/7 go 
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The latter can be simplified by noting 

go(0) = (Ao - Bo)"'j(Ao + Bo), (70) 

with 

Therefore, 

= (Ao - e^2<^%o)"'j(Ao + e^^^^Bo), (72) 
and defining, by analogy with fITT]) . 



A.-e-B„^(-|t y, A„.e--B„^(^-|* j) . (73) 

then the matrix g^f^ can be expressed in exactly the same structural form as go(0) in 
terms of 3 x 3 matrices, as 

(0) ^ ( ^ -Q^'?* ^ Z'^j' \ (74) 
l^.(T,-R,+Q-R,) RjQ-j' ^^'^ 

where the tt— periodic submatrices are 

Q<f> = = cos^ Q + sin^ 0T + i sin0cos0(R — R"^) = 

= T+ = cos^ T + sin^ Q - 2 sin cos (R - R+) = (75) 



R0 = cos^ R + sin^ R+ + zsin0cos0(Q - T) = R^ 



+ 



The submatrices defined in eq. (!65l) therefore have alternative integral expressions, from 
dSHD and dZH]), 

S = - (Q^^R^), H = - (Q-i> , B = (T^ - RJQ^'R^). 

Crucially, Q,^ is positive definite for n > 1. In order to see this, first note the obvious 
positive definiteness > if sin0 = 0. For sin 7^ we have 

= -sin2 0A( -icot0) where A(A) = A^Qo + A(Ro - R^t) - To, (76) 

which is positive definite by virtue of the fact that det A (A) = does not admit pure 
imaginary roots for A [13, §3.2.1]. By the above arguments, the matrices H and B are 
negative and positive definite, respectively. Consequently, H and B are invertible, and so 
the identity (!66|) 9 implies that the matrices I — and hence I — S"*"^ are also invertible. 
Note that the positive definiteness of B confirms that det B 7^ for n > 1, which is when 
there is no stress- free modes. The cases n = 0, 1 are discussed in ^5.2[ 
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5.1.2 The impedance matrices Zq and Zonp 



We are now in a position to express the impedance matrices in terms of the integrals. As 
before, we set a = 1, 2, 3 and a = 4, 5, 6 for the physical and nonphysical triplets, respec- 
tively. Inserting Li = — zZqAi, L2 = — zZonpA2 in fHSjl o and using the same argument as 
in [l| to maintain that det Ai 2 7^ implies 



I I W Ai 

-iZo -iZonp) V A2 

which, together with eqs. (IMj) (|^ . yields the matrix identity 

S iH \ / I I \ _ / I -I 

iB -S+y l-iZo -iZonp/ V~^^0 «Zonp 



(77) 



The first line yields explicit expressions for the impedances Zq and Zonp, 

Zo = - H-^S, Zo,p = -H-i - H-^S, (79) 
and the second line gives the equivalent expressions 

Zo = - (I + S+)"'b = -B (I + S)'\ Zo,p = (I - S+)~' B = B (I - S)-^ . (80) 
Hence, Zq and Zonp are Hermitian by virtue of (!79|) and (!66|) . 

Zo = Zq , Zonp = Zg^p. (81) 

Using ( !79|l and (165|1 implies 



Zq — Zonp — 2H 



-i = z(AiA+)-'(=-z(A2A+)-') (82) 



in agreement with the formula in §4.41 

The representation (!79|) of the physical and non-physical central impedances is similar 
to that of, respectively, the non-physical and physical half-space impedance Z [v) in [l|; 
however, the two addends in the right members of ( !79|) are generally not the real and 
imaginary parts of the impedance as is the case in the expressions of [l|. Consequently, 
the non-physical central impedance is not minus transpose of the physical one, unlike the 
half-space impedance where Z (f ) = —7j^^{v) [H. Another noteworthy difference is that 
the expressions similar to ( IHOj) are not unreservedly valid for the dynamical Z {v) because 
the analogue of B has zero determinant at the Rayleigh speed. 

The above results can in the main be obtained by following the alternative method of 
deriving the Lothe-Barnett integral formalism which was proposed by Mielke and Fu in 
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5.2 Remarks on n = 0, 1 

Let us now discuss the implication of the special cases n = 0, 1 in the present context. 
Occurrence of a non-semisimple go (0) due to a pair of degenerate eigenvalues A'-^'^-' = 0, 
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which are spht between physical and nonphysical triplets, makes the cases n = 0, 1 
tantamount to the limiting state v = v oi the elastodynamic problem for a half-space 
[H, HI, mi (more specifically, to the so-called exceptional limiting state in view of the zero- 
traction mode corresponding to A*^*^'^) = 0). The Lothe-Barnett integral formalism on the 
whole is well-defined ai v < v. Any difficulties occurring a.t v = v are due to the non- 
integrable divergence acquired at f = w by the angularly varying Stroh matrix N [v, (p) , 

which is a counterpart of glf^ ■ A similar exception arises with g^f'' for n = 0, 1 due to 
Q^^. The argument underlying positive definiteness of for n > 1 no longer applies 
for the cases n = 0, 1 which admit rigid-body motion. This can be ascribed to the fact 
that det k^^-^) = in (n^ - 1) = so that det A (A) = has the root A^"-^) = 0, where A (A) 

is given by (f76|) . Thus Q,^ for n = 0, 1 is positive semi-definite, with det Q,^ = at 
coscj) = 0. That is why the cylinder's version of the integral formalism cannot generally 
be extended to the cases n = 0, 1. The exception when this is yet possible is the case 
n = for a cylindrically monoclinic material with the symmetry plane orthogonal to the 
2;-axis. This case simplifies due to the simultaneous occurrence of as the null vector 
of K*-*^-* and of the uncoupling of the 2;2;-components. Hence the upper 2x2 blocks of the 
integral-formalism relations remain valid. Such a state of affairs also has a direct analogy 
with the theory of surface impedance in half-space, namely, with the case of a symmetrical 
sagittal plane, which is when the in-plane modes are unaffected by the limiting state vsh 
of the uncoupled shear-horizontal mode [13, [HI . A careful remark is in order regarding 
Eq. (I80ll . For n = 0, the rigid-body displacement corresponding to A*^"-* = 1 and parallel 
to eg is the null vector of B and the eigenvector of S with the eigenvalue (A*-*^-*) = 1. Hence 
the upper 2x2 blocks of B and I — S, I — S"^ are singular, whereas those of I -|- S, I -|- S"*" 
are not. Thus flHOj) is not valid for Zonp even for the monoclinic n = case. 

Finally, it needs to be added that analysis of the half-space integral formalism asv 
[11, m, [29[[39| shows that although the formalism diverges on the whole in this limit, the 



integral expressions for the surface impedance Z (v) remain well-defined at v = v. This 
asymptotic property of Z (v) is not however directly relevant to the central impedance Zq 
of a cylinder, since the diverging cases n = 0, 1 cannot be approached 'continuously in 
n\ A more appropriate treatment is either asymptotic analysis of Z (r) — )■ Zq as r — )■ 
or else other, explicit, methods of deriving Zq for n = 0, 1, see §5.41 



5.3 Definiteness of Zq for n > 1 and semi-definiteness for n = 0, 1 

It has been noted above that the structure of the physical and non-physical central 
impedances ([79l) resembles that of, respectively, the non-physical and physical surface 
impedances [l| for a half-space. This suggests the inverse correspondence of their sign 
properties. We will outline a formal proof. Similarly to the Lothe-Barnett theory for 
a surface impedance, the insight does not follow from the integral formalism but relies 
instead on static energy considerations. 

As before, assume first that n > 1. The central impedance is related to go (0) inde- 
pendent of uj and k^, therefore we can invoke the 2D static solution A/" (r) = Fdiag (^r'^"^ 
(cu = 0, kz = is tacit below). The time averaged energy associated with these solutions 
is [10, Eq. (21)] 

= - — ^(U+V-V+U). (83) 
8r dr ^ ^ ^ ^ 
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Inserting the physical solutions with eigenvalues ReXa > and eigenvectors 7oq = 
(a„, \a)^ (a = 1, 2, 3) of go (0) from §4.1.11 and using the central impedance Zq = Zq 
leads to 

Wrdr = -^Yl (^2^°^" ~ a;Zoa« > 0, Vrs > ri. (84) 

Hence Zq for n > 1 is negative definite. The same consideration using the non-physical 
solutions with ReA^ < 0, a = 4, 5, 6, implies that Zonp for > 1 is positive definite. As 
expected, this is opposite to the sign properties of the physical and non-physical surface 
impedances Z (u), Z^p (v) for the static limit v = 0. 

In the case n = 1, the above proof applies unchanged for the physical Z^^ except that 
it is negative semi-definite due to the presence of the rigid-body motion mode. In the 
case n = 0, the same conclusion of negative semi-definite Zq°^ follows from an explicit 

calculation of Zq^'' = limr-->o Z^^^ (r) presented in §5.4.21 It is noted that the rigid-body 
displacements causing 

det Zo = for n = 0, 1, (85) 

are related to the existence of low frequency (long wavelength) guided waves in rods: 
longitudinal, torsional (n = 0) and flexural (ra = 1), see §4.1.21 and j20| . 

5.4 Explicit expressions for the central impedance matrix 

In this subsection we develop other procedures for determining Zq, including for the special 
cases n = and n = 1. 



5.4.1 Zo for n > 

The central impedance Zo for n > is defined by any of the optional relations (including 
( IS^ t) that may be written similarly to f l53|) as 

L = -zZoa„ (a = 1,2, 3) ^ 1 = -zZoa ^ Zq = iLiA^\ (86) 

T T 

where 7 = (a, 1) is an arbitrary superposition of the physical eigenvectors 7^ = (a^, !„) 
of go (0) with a = 1,2,3. The matrices Ai = ||ai,a2,a3|| and Li = ||li,l2,l3|| may be 
related to one another using identities such as (24) and (26) of |20|, 

L = i(A„Qo + Ro)ac, (a = 1,2,3), 

implying 

Li = z(QoAiA + RqAi), where A = diag(Ai, A2, A3). (87) 

The eigenvectors a^ (a = 1,2,3) are null vectors of A(Aa), see. eq. fl76|) . and conse- 
quently, 

QoAiA' + (Ro - R^)AiA - ToAi = 0. (88) 

Equations ( l87|) and ( !88l) provide a pair of expressions for the central impedance, each in 
terms of the displacement eigenvector matrix only, 

Zo = -Ro - QoAiAA^i = -R+ - ToAiA"^A^\ (89) 

The freedom afforded by these simultaneous identities will prove to be useful when ma- 
terial symmetry reduces the matrix size to 2 x 2, see §7.1.11 
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5.4.2 The central impedance for n = 

The algebraic Riccati equation f lSI?]) for = leads to a constructive solution for Z^^\ 
which must satisfy 




r (Cl2 -Cl6 0\ -X /-C22 C26 0\ 

Qo^<^ Z(°)+ C26 -C66 \+[ C26 -Cgg =0. 

\C25 -C56 0/ J V 0/ 



Noting that 




\ _ /O ci6 0^ 

Cl6 C66 C56 Qo ^ = Qq ^ ^66 
/ \0 C56 0, 

it is clear that the solution of the algebraic Riccati equation is of the form 

Z{,°) = I I , (90) 




where the scalar z^^^ satisfies a quadratic equation 

(^("^e^ + p)'^Qo^(2;(°)e^ + p) - C22 = 0, with = [c^, Cae, C25). 

The physical root must have negative real part in order to be consistent with eq. (I92 
below, implying 



z 



(0) 



(Qo')n 



[ - - V + (c22 - P^q) {Qo ^) iJ ' where q = Qq ^P- (91) 



An alternative method is to take the limit r — )■ of the solution for jV of |20|, Eq. 
A4]. The result is again eq. ( I90l) where, by definition, 2;^°^ is (i^r) times the ratio of radial 
components of the traction and displacement of the eigenvector 7*^°-' of go(0) corresponding 
to its eigenvalue X = k. These were found by Ting 2J], from which z^'^^ follows as 



.<») = -iI±^, (92) 

C55C66 ~" 



where Q = det Qo and, using the notation of [24 

C12 C26 C25\ / C22 C26 C25' 

= det I C16 C66 C56 , Y = det I C26 cee cse 

^Ci5 C56 C55/ yC25 C56 C55^ 

Equivalence of the expressions f l9T]) and fl92|) follow from identities such as qi = W/ Q and 
(C22 - P^q)(Qo')n = (^^ - W^)/Q\ Note that QY > [H, Eq. (B2)]. 
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5.4.3 Zo for n = 1 

The physical triplet of eigenvalues Xa^ and eigenvectors 7a'' = \a^)'^ of go (0) (a = 
1, 2, 3) includes A'-^-' = 0. It corresponds to a rigid-body rotation about the ^— axis with 
displacement vector a'-'^^ = (l,z,0)"'" and zero traction l^^^^ = 0, see (lol, Eq. (52)]. Hence 
by dSlDa Z[,^^a(i) = 0, i.e. a*^^) is the null vector of Zg^-* for any anisotropy. This property, 
combined with (185|) . implies that Zq^^ has the structure 

Zq^"* = I — m a —ic I with I I negative definite. (93) 




The 2x2 matrix becomes diagonal (c = 0) for symmetry as low as monoclinic and an 

-'0 



explicit form of Zq^^ can then be found, see §7.1.11 



5.5 The matrix Zq at large azimuthal order n 

For n ^ 1 we assume an asymptotic expansion of the impedance in inverse powers of n: 

Zo = nzo + zi + n~^Z2 + . . . , (94) 

where zo, zi, . . . , are independent of n. Substituting into flS^ and comparing terms of 
like powers in n yields a sequence of matrix equations, the first of which is 

(zo - iR^) Qo ' (zo + ^Ro) - To = 0. (95) 

This algebraic Riccati equation can be identified as eq. flTTl) with system matrix Q = zfcgN 
where kg = n/r and N is the (static) Stroh matrix for the sagittal plane defined by 
e^, eg = n, m. The subsequent identities are inhomogeneous Lyapunov equations 

E+z, + z,E + f,(zo,zi,...z,_i) =0, J = 1,2,..., (96) 

with the constant matrix operator 

E = Qoi(zo + zRo), (97) 

where 

fi = zToK + zKTo + E+RoK + KRj^E, etc. (98) 

The leading order impedance zo is the solution of the matrix algebraic Riccati equation 
fl95l) . and may be determined by the methods discussed above (via the eigenvectors & 
eigenvalues, or the integral representation). Subsequent terms zj, j = 1,2, . . ., satisfy a 
Lyapunov equation fl96|) with different right hand sides but the matrix Lyapunov operator 
is the same in each case. The solution of the Lyapunov equation depends upon the 
spectrum of E, and since the eigenvalues of E have negative real part, it follows that the 
unique solution is 

oo 

zj = J dse'^^i.e'^. (99) 



The asymptotic sequence in inverse powers of n can thus be evaluated to any desired 
order. 
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6 Radiation impedance matrix 



The radiation impedance is relevant, for instance, in a configuration of infinite outer extent 
in which the cyhnder is inhomogeneous in r < tq for finite tq and uniform otherwise. It is 
always possible to split the linear total field into incident and scattered components, such 
that the scattered solution in r > tq has only positive radial energy flux. The radiation 
impedance is defined by the subset of wave solutions with this radiation property. 

6.1 Explicit form of Zrad(^) 

An alternative partitioning of the integral matrix is required to account for the separate 
radiating, i.e. non-zero flux, modes. This may be accomplished by a change of basis that 
brings about the diagonal form of the flux matrix jV"*" (r) TAf (r) [l2|. Proceeding from 
the integral matrix J\f which satisfies ( l47l) . and hence composed of modes with zero radial 
flux, we first convert T to diagonal form by an orthogonal transformation, 

T = WJW+, J=(^-^ J), W=i=(_^j (W+W = I). (100) 

Then referring to the notation fH9l) i for A/" satisfying eq. fHTl) . we are led to the following 
partitioning of the integral matrix, 

/U+ U_\ /U1-U2 Ui + U2\ 

The + and — suffices indicate modes that have positive and negative flux in the radial 
direction, which is evident from the sign of the flux defined by eq. (!27|) . and the flux 
condition (1471) which becomes 

Aft (r) TAT I (r) = J forn > 1 Afi (r) JAff (r) = T) . (102) 

Extension of this identity to the special cases n = 0, 1, is contingent on the details of the 
Frobenius solutions, see §4.1.21 In the cases of transversely isotropy and isotropy, the + 
and — modes correspond to radiating (outgoing) and incoming Hankel function solutions, 

H^'^ and Hi'^ respectively. 

The wave-based partition (llOip provides the required modes to express the radiation 

impedance, defined in ( 139|) . with Uscat, Vgcat — ^ U+, V+, 

Zrad(r) = zV+U;\ r > 0. (103) 
It is important to note that the radiation impedance is not Hermitian, since according to 

(TO 

Zrad - = -Z ^ 0, (104) 

which implies in fact that i(Zrad — Z^^) is Hermitian and positive definite. 

As an example, consider SH wave motion in a uniform isotropic solid with kz = 0, 
for which the scalar radiation impedance is Zj.s,d{r) = —c^^krHn^ {kr)/Hn\kr) where 



28 



k = a/p/ C44 and Hn '^ is the Hankel function of the first kind, see eq. (I119p . Using 
known properties of cyhndrical functions yields for this case 

z(Zrad - Z+J = An-'c,,\Hi'\kr)\-' > 0. 

Note that for large values of kr the SH radiation impedance is ^radl'") = —ikrc44, + 
ic44+0((A:r)-i). 



6.2 Asymptotic form of Ziad(^) as r — > 00 

Assume that the cylinder material is homogeneous for tq < r < 00, for some finite radius 
Tq. As r — 00 the impedance Zj-adl'"); which we recall is defined with generalized traction 
vector V = irT(r), may grow without bound while r~^Zi.ad('") tends to a planar limit. 
This behaviour is evident for the SH radiation impedance considered in §6.11 which is 



proportional to r as r — )■ 00. We therefore assume that the radiation impedance has the 
form 



'Z'rs.dir) = k^rZoo + 0(1), r cx), 



(105) 



where k^Z^o is a constant matrix. This can be found by considering the large r limit of 
the differential system ( 1221) . which reduces to its plane wave asymptote with r playing the 
role of a rectangular coordinate: 



— 0(r) = ifo0(r) 
dr 



where ^(r) is a 6— vector and fo a matrix constant [12 
<^W=(v(rp' y(r)=ik:'r(r), fn = fc. 




(106) 



(107) 



The six independent solutions to fll07p may be separated into triplets according to their 
flux properties, with U+, V+ signifying the outgoing, or radiating solutions. The limiting 
radiation impedance is then defined by analogy with fllOSp as 

Zoo = ^V+U;^ (108) 

Properties of Zoo can be deduced by noting that the system f ll06p is equivalent to that 
for a half-space with the identification fo = A;^N(w), where v = u/kz and N (w) is the 
elastodynamic Stroh matrix for the sagittal plane {e^, e^} = {n, m}. This enables us to 
equate the limiting radiation matrix Zqo with the surface impedance matrix Z {v) for a 
homogeneous half-space 0]. Consequently, Zoo = "^oo subsonic v, i.e. < v < v. The 
possibility of Zoo being Hermitian seems at odds with the conclusion f ll04p : however, it 
should be borne in mind that Zoo is only the leading order term in the asymptotic series 
implicit in f llOSp . The subsonic situation may be understood in the context of the SH 
radiation impedance example above with the wavenumber k formally taken as imaginary, 
in which case the Hankel function is replaced with the modified Bessel function of the 
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second kind via the identity Hn \x) = 27r ^{—i)'^'^^Kn{—ix). Conversely, Z^o is not 
Hermitian for > {) jlj. The equivalence with the half-space problem also implies that 
Zoo is a solution of the algebraic matrix Riccati equation 

(Zoo - ^Pf) Qc' (Zoo + ^Pc) - M, + p,vH = 0, (109) 

where the suffix c indicates the constant values in r > tq. Equation (11091) can be deduced 
by analogy with eq. (p5|) . noting the presence of the additional dynamic term pcv"^! in fg 
and hence in (I109p . The Riccati equation indicates that as /c^ — ?■ the matrix hzZ^o 
Zooo where ZoooQ^^Zqoo = — Pci^^I, with a unique solution satisfying (I104p . and hence 

lim r-^Z,ad(r) = -lupl^^Ql^^ for = 0. 

r— >oo 

Note that taking Qc with cis, csg = and C44 = C55 factors out the asymptotic form 
Zoo = —ikrc44 of the above-mentioned scalar radiation impedance Zy-^d for the SH waves 
in an isotropic solid. 

Finally, it is emphasized that developments in this subsection are irrelevant to the 
solid- cylinder impedance Z (r) which, by construction, is Hermitian at any r and for any 
f (= uj/kz). This in fact implies that r~^Z (r) cannot become constant as r — 00 because 
otherwise the arguments subsequent to eq. (11051) would violate the unconditional Hermitic- 
ity of Z (r). For instance, the out-of-plane impedance, Z^ir, 0) = — C44 k2rJ'^{k2r) / Jn{k2r) , 
see eq. ([US]) 2, has no large-r limit. 



7 Explicit examples of the solid impedance 

The central impedance Zq is first presented for several cases of material symmetry, includ- 
ing monoclinic and orthorhombic. A semi-explicit form for Z(r) is possible if the material 
is transversely isotropic, providing a check on the numerical calculations in §7.31 



7.1 The central impedance Zq 

It follows from its definition through go(0) that Zq depends at most on 15 of the 21 
possible elastic moduli. The six redundant moduli are those with suffix 3 occurring in the 
Voigt notation. 



7.1.1 Monoclinic symmetry 

For monoclinic symmetry with the symmetry plane orthogonal to the 2— axis the impedance 
has the structure 

7 ° 



;iio) 







where Z±q and Z^o are the in-plane and out-of-plane impedances, respectively. The out- 
of-plane scalar impedance follows from [20|, Eqs. (37), (38)] as 



^20 



nJcMC55-cl^- (111) 
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For n > 1 the in-plane impedance can be expressed in semi-explicit form in terms of 
the eigenvalues Xj, ReXj > 0, j = 1, 2, of the 2x2 matrix g_Lo(0) formed from the upper 
left block of go(0). By use of the following identity for 2x2 matrices, 

A + = (Ai + Aa)! (AiAs^O), 

the formulae in (!89l) may be combined to eliminate the explicit dependence on the eigen- 
vector matrix, with the result 

Zxo = -^(Ro + R^) - (Qo ' + A1A2T0-1)-' [(Qo-i - AiAsTq-^) ^(Ro - R+) + (Ai + X^)!] . 

Note that the matrices on the right hand side are all 2x2, i.e. Rq = R±0; etc., and 
the eigenvalues are the two roots of the quartic det A_l = from eq. ( ff6l) with positive 
real parts. The block impedance Z^o depends upon the six in-plane moduli, c^s (/^, <^ = 
1,2,6). 

For n = 1 the in-plane impedance possesses a null vector as described in §5.4.3[ and 
based on the required Hermiticity, it must have the form 

ZS = .«(_\ ;)=.«e+e, e=(l,z). (112) 

The algebraic Riccati equation (159|) then reduces to 

{(^We-zeR;^)Qo^(^We+ + zRoe+) - eToe+}e+e = 0, 
implying a quadratic equation for z^^\ 

(cii -|- Cqq) — 2z^^^ (ciiCge — c^6 ~ C12C66 + C16C26) 

- (C11C22C66 + 2Ci2Ci6C26 - C11C26 - C22C^6 ~ ^\2^&^ = 0- 

The unique physical ^(1) is, according to §5.3[ provided by the negative root. We note that 
the eigenvalues for the in-plane modes are A^^^ = and X^2^ that is the physical (positive 
real part) root of 

A^ (ciiCge — c\g)+2iX (ciiC26 — Ci2Ciq) + {ciq — C2q)'^+Cqq (2ci2 — Cii — C22)— CiiC22+c^2 ~ 0- 

(113) 

For n = the general expression f PT]) reduces to 

C66 y V Cge/ V ^66/ 

7.1.2 Orthorhombic and tetragonal symmetry 

For the orthorhombic symmetry and n > 1, the in-plane impedance z||"q is given by the 
upper 2x2 block of (EH]): 

Z.o = ( - f^- ) A, ('^ ) n > 1, (114) 

\-tnceQ Cqq J \0 ceej \0 X2J ' 
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where Ai 2 are the physical roots of the equation 

A^ciicee - [n^ (C11C22 - 0^2 - ^cuCae) + cqq (cn + C22)] + C22C66 [n^ - l)^ = 0, 

and A_L = ||ai^,a2_L|| is composed of the null vectors of the matrix A_l (A), and can be 
expressed 

A?C66 - C66 - n'^C22 -m[A2(ci2 + Cqq) - C22 " Cqq] 

-m[Ai(ci2 + cee) + C22 + cqq] Xjcn - C22 - n'^cm 

For n = 1 the scalar in-plane impedance is 

^(^)-^^(cn-c,2-CnAW), 

Cll + C66 V / 



where A^^^ = V {cnC22 - ^2 + cnCge + C66C22 - 2ci2C66) /(cnCge) is the (physical) root of 
(I113p simplified for the orthorhombic case. 

For tetragonal symmetry with Cig = C26 = the in-plane impedance is unchanged 
from f llMp . and the out-of-plane impedance f illip further simplifies due to C44 = C55 (on 
top of the orthorhombic condition C45 = 0). 

7.1.3 Transverse isotropy and isotropy 

The central-impedance matrix reduces for transversely isotropic symmetry to 





1 2C66(-™) 
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—ncM/ 





(115) 

'1 0^ 



ZP = -2(cii - C66) I 

.0 0> 



which applies, of course, to isotropy (C44 = Cqq). Equation (11 15p is also derived in the 
next subsection. 

7.2 The solid-cylinder impedance Z(r) for transverse isotropy 
7.2.1 General formulation 

The constitutive relation for combined with eqs. fl2T]) and ([53]) 3, implies for any material 
anisotropy, 

Z(r) = -R-zA;^rP-Q(r-^Ui)U7\ r > 0, (116) 

dr 
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where the matrix Ui (r) is any unnormahzed triad of independent physical solutions. The 
radiation impedance is obtained if the matrix is replaced with U+(r) comprising linearly 
independent radiating solutions. The difficulty in applying flll6p is that explicit matrix 

solutions for Ui or U+ are not generally available except under certain restrictions on 
material symmetry, such as transverse isotropy. 

Assuming transverse isotropy, solutions for the displacements that are either regular at 
r = or radiating to infinity can be constructed in terms of cylinder function by adopting 
the representation of Buchwald [40| (see also [HI). Thus, 

-Cn{k3r)\ 



V(r) 



m 

k^r ^ 



m 







where the principal wavenumbers ki, k2, ks, and auxiliary wavenumbers Ki, K2, are 
y/a"^ — b ,2 pu'^ — c^^kl CQek^ — Cukf 



kf 



2CiiC44 

(Cii + C44)pu;^ + 



kl 



^13 



C66 

2C13C44 — CiiC33)/c 



kz{Ci3 + C44) 

b = AciiCii{pijj'^ 
(1) 



1,2 



C33kl){pu'^ - c^kl 



and Cn = Jri for displacements regular at r = , C„ = Hn for radiating solutions, where 
Jn are Bessel functions and Hn^ are Hankel functions of the first kind. 

Evaluating (11161) and simplifying terms using the identities C44K1K2 + c^^kl = 0, 
cii(ki/c| — K,2kl) = cee^K^i — ^2), we find that the solid cylinder impedance for transverse 
isotropy is 



Zir) 
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with the non-dimensional quantities 
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(117) 



6(62/1 - 62/2) - n'^iVi - 2/2) 



kivSM (J = 1,2,3). 



Cn{kjr) 

The central-impedance limit may be extracted from f lll7p by writing it in block form 

Z±(r, kz 
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where the dependence on both r and is emphasized. For /c^ = we have kj = u/cj 
with pcf = Cii, pc^ = C44, pel = Cqq, and the impedance reduces to 



Z(r,0) 



Z^(r, 0) = 2q 



'66 



Z±(r,0) 


1 in} 

-in 1 






+ C66(/C3'^)^ 



with Zz{r, 0) = — C44 A;2'" 



(119) 



'^3' C„(fe3r-)/ 



Taking the hmit r — )■ of (11191) with the interior cyhnder functions C„ = Jn gives eq. 



7.3 Numerical example 

A procedure was outhned in §4.31 for calculating the solid-cylinder impedance using two 
separate numerical solutions. The Riccati equation ( 158|) is first integrated starting from 
r = with the central impedance matrix Zq as initial condition. The integration proceeds 
up to r = ri where ri lies below the first singularity of Z(r). For r > ri the impedance 
is obtained from ( MT\i o as the solution of the the matricant-based system ( 160|) . with the 
Riccati solution at ri serving as the initial condition. As an illustration of its practicality, 
the two-stage algorithm was implemented with representative results plotted in Figure |2l 

The initial step in the computation requires the value of the central impedance, which 
was calculated using eqs. ( I90l) and (19T1) for n = and the formula Zq = zLiAj"^ for n > 0, 
see (15^ 9. with Ai, Li defined by the numerical spectral decomposition of go(0) and the 
appropriate selection of its three eigenvalues with positive real part. It was confirmed 
that the computed Zq satisfied the algebraic Riccati equation (|59l) . with error always less 
than 10~^^. Numerical integration of equations ( l58l) and (|60l) was accomplished using the 
Runge-Kutta (4,5) routine in Matlab. In order to assess the accuracy of the numerical 
results the computed matrix Zcompi'f') and the analytical solution for Z(r) of (11171) were 
compared. For the examples shown in Figure [2] it was found that the spectral norm of 
the difference satisfied ||Zcomp('^) — Z(r)||2 < 10""^ at all points. The curves in Figure [2] 
use r = 1 as the 'cross-over' coordinate, but similar accuracy was found for other values 
as long as they lie below the first singularity of Z(r), which for the parameters considered 
is > 2. In all cases the transition from the Riccati to the matricant based solution was 
found to be smooth. 

This numerical procedure is designed to handle the coordinate-based singularity present 
in the system equations ( l22l) at r = 0, and can be continued, in principle, to any finite 
r. At the same time the computed impedance Z(r) will grow without bound at discrete 
values of r > associated with waveguide modes of the traction-free cylinder. The point 
of the algorithm is that it will continue to provide accurate solution regardless of the 
presence of two distinct types of singularity at r = and at finite values. 
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Figure 2: The curves show {n^ + 1)^^| det Z(r)| for n = 0,1,..., 5. The material is 
isotropic with {cii,C66,p} = {4,1,1} and {uj,kz} = {1,0.2}. The Riccati equation flSSl) 
for Z(r) was integrated to obtain the curves for < r < 1, starting from r = with the 
known Zq of eq. f lllSp . For r > 1 the system equations ( 125|) were integrated and eqs. ( 160|) 
and ( l6Tl) used to find Z(r), starting from the Riccati solution at r = 1. 



8 Conclusion 

Impedance matrices appropriate to cylindrically anisotropic radially inhomogeneous elas- 
tic materials have been defined and procedures for their determination developed. In 
the process a new impedance matrix has been revealed as of central importance for wave 
motion in cylinders with on-axis material. The solid-cylinder impedance matrix is a char- 
acteristic property of the cylinder, with no free parameters apart from frequency and axial 
wavenumber. The impedance may be defined as the unique continuation of its on-axis 
limit, the central-impedance matrix, which is a simpler object dependent only on (a sub- 
set of) the elastic moduli. Two methods have been described for constructing the solid 
cylinder impedance at r > 0, one based on a Frobenius series solution, the other using a 
differential Riccati equation. In addition to providing practical means for computation, as 
has been demonstrated for the latter approach, the methods shed light on the structural 
properties of the impedances. The Frobenius solution offers direct proof of uniqueness 
and Hermiticity, while the Riccati solution provides a stable method to integrate the 
otherwise singular system of equations at r = 0. The radiation impedance matrix, suit- 
able for infinite radial domains, has been defined and its properties delineated. We have 
found it instructive to compare the cylindrical impedance matrices with the surface wave 
impedance for a homogeneous half-space. The central-impedance matrix is the negative 
semi-definite counterpart of the static surface impedance, and the large r limit of the 
radiation impedance is closely related to the surface wave impedance with v = u/kz- 
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One purpose in developing these impedance matrices is the significant advantage of- 
fered by the impedance approach in solving boundary value problems. The solid-cylinder 
impedance matrix provides perhaps the simplest method to arrive at the dispersion equa- 
tion of a radially inhomogeneous solid cylinder. In this regard we note that, by analogy 
with the conditional (3 x 3) and two-point (6 x 6) impedances of an annulus [l2|, the eigen- 
values of the solid-cylinder impedance should be monotonic in u: at any fixed kz-, which can 
be helpful for finding numerical solutions of the dispersion equation. In a wider context, 
the impedance matrix in conjunction with the radiation impedance, can serve in formu- 
lating scattering of acoustic and elastic waves from solid cylinders. Other applications 
that we envisage include the use of impedance matrices for solving problems with dis- 
tributed forces within the cylinder, and applications involving 2D-inhomogeneous or lat- 
erally bounded planar and cylindrical waveguides jH], |43j where the algebraic impedance 
matrices discussed here become differential operators. 

Another no less important reason for investigating the impedance matrix in the cylin- 
drical context is that it affords new insights on the nature of elastodynamic solutions in 
anisotropic elasticity. It is remarkable, for instance, to find the Riccati equation appear 
as a natural method for solution in cylindrical elastodynamics. The Riccati equation, 
in fact, implies that the central- impedance solves an algebraic Riccati equation, which 
in turn leads to direct methods for its evaluation using analogies with the surface wave 
impedance. Differential Riccati equations have been found useful in few elastic wave set- 
tings, e.g. ssimiiii. Its appearance here suggests it has wider potential application 
in computational elastodynamics. 
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